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ABSTRACT 


Two approximate parametric interval estimation methods for system rehabilitv using 
component test data are developed and evaluated. One method can be apphed to any 
coherent system with components which have exponential failure times with possibly 
different failure rates and different mission operating times. This method estimates the 
ratios of component failure rates which are then used to develop the approximate lower 
confidence limit. These ratio estimates are developed with and without jackmife methods 
and the two results are compared. This procedure is verv accurate and simple to com- 
pute, requiring the use of standard chi-square tables. This razio method is subsequently 
extended to coherent svstems with components whose failure times have a HV eibull dis- 
tribution. A nearly exact parametric lower confidence lumt for PLY > x) is developed and 
evaluated where x 1s given and X has a normal distribution with unknown mean and 
variance. This procedure is also simple to evaluate and requires the use of Student t ta- 
bles. 
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I. BACKGROUND 


Parametric confidence interval procedures for the reliability of mechanical systems 
are much less developed than procedures for electrical svstems. This is due to the more 
complicated failure distributions used to model mechanical hardware. The failure rate 
of a series system of independent components, with exponentially distributed failure 
times and equal mission times, is the sum of the failure component failure rates. This 
property has permitted development of numerous methods for system reliability of series 
and other coherent systems using component failure data. The Heibull and extreme 
value distributions have been used in life testing methods for both electrical and me- 
chanical devices. Several methods have been used for obtaining point estimates of pa- 
rameters for these distributions and for the reliability function itself, Harter and \loore 
Beet. 1. pp. 889-901], Mann [Ref. 2: pp. 231-256]. 


The derivation of simple confidence limits for the rehabilitv function for the extreme 
value distribution with parameters z and 0 has posed problems. because methods based 
on a pivotal quantity such as (z — z)/6 are inadequate. Lawless [Ref. 3: pp. 355-3064]. 
Johns and Lieberman [Ref. 4: pp. 135-175] and Thoman. Bain and Antle [Ref. 5: pp. 
363-372] have developed nearly exact procedures for confidence limits for the rehabilty 
functuon of the Weibull and extreme value distributions. Schneider and Weissfeld [Ref- 
6: pp. 179-186] have developed interval estimation methods for percentiles of the TH esbull 
and extreme value distribution based on censored data. Although extensive methods 
have been developed for interval estimates of the rehability of a single component with 
the Heibull or extreme value failure distribution, very few parametric interval methods 
have been developed for system reliability using component test data with JVeibull failure 
distributions. 

Two approximate parametric interval estimation methods for reliability of coherent 
systems using component test data are developed and evaluated in this thesis. Evalu- 
ations of these procedures are performed using computer simulation for series svstems 
onlv. One of these methods is developed for the reliability of a series system whose 
components have exponential failure distributions and different mission times. This 
procedure was found to be quite accurate and can be applied to coherent systems in 


general. This first procedure is then extended to the case where components of the svs- 


tem have HWerbudl failure distributions. The method used is an extension of a non- 
parametric method developed by Myhre, J., Rosenfeld, A. and Saunders, S. [Ref. 7: 
pDeanoee 2 

The normal distribution is used extensively in some mechanical reliabilty models. 
Maximum likehhood and minimum variance unbiased estimators for PLY >.x,) when 
both uw and o? are unknown were developed over thirty years ago, Lieberman and 
Resnikoff [Ref. 8: pp. 457-516], Folks and others [Ref. 9: pp. 43-50]. and Barton [Ref. 
10: pp. 227-229]. Exact interval estimation procedures for P(X > x,) were developed by 
Owen and Hua [Ref. 1]: pp.285-311] using the non-central ¢ distribution. Letting X de- 
note strength and Y denote stress, mechanical reliability is sometimes modeled as the 
value for PLY > ¥). Approximate interval estimation procedures for this probability 
when Y and Y are assumed to be normal have been developed bv Church and Harris 
(Ref. 12: pp. 49-54] when the mean and variance of } are known. Downton [Ref. 13: pp. 
5581-558] modifies their procedure slightly to get more accurate bounds and suggests an 
approximate procedure when the means and variances of both VY and ¥ are unknown. 
Lower confidence intervals for PLY > ¥) obtained under the assumption of normality for 
A and Y can lead to serious error when either .Y or ¥ or both are truncated well into the 
tails. Consequently, P(X > x,) may be a more reasonable model of mechanical reliability 
where x, 1s chosen conservatively. A very accurate approximate lower confidences 
procedure for P(X > x,) 1s developed and evaluated in this thesis. It can be computed 


easily. 


If, INTERVAL ESTIMATION PROCEDURE - EXPONENTIAL CASE 


A system of independent components is coherent if an increase in reliability of any 
one of its components does not cause a degredation in system reliability. Suppose a 
coherent system has & components. We assume that the failure distribution of compo- 
nent / 1S exponential with failure rate 4,. Then system reliability R, can be written as a 


fomeon Of z,¢%, f=1,2,°**,k; iLe., 
R(t) = gls;, Ad, errs | cos i, (5, eee, ,) 


where ¢ = 2,1) is the operating time for component /, / =1, 2.e+¢+,k., 


wet — Max 73,8, /,, andy =///A, (= 1) 2.<<e«, k. When one can write 


Ks(1) a gua ESO RA se ty) 
fame y .were known and /,,,,,, were an upper 100(1 —a@)°> confidence limit for /,, , the 
corresponding lower confidence limit for R(t) would be 


A 


ROL) a ACs UA Oe atin. 95 ty) 


Specifically, if the svstem is a series system of independent components, so that 


k n 
R(t) = ost _ » ‘tie = exp, — Aa rt, (2.1) 
i=) i=] 


then, 


fn 

A A 

Rs(Qpa) = exp, — ere Vij 
i=! 


If 2, items of component / are tested until failure, 7, denotes the total test tume accumu- 
b b 

Fred by all vx. items andv=>n7,then 27,577, 1s 73,. See Bain and Engelhardt (Ref. 
= (=1 

14). 


An upper confidence limit for «,, 1s 


2 
A Xx, 2n 


Am, U(2) k 
2 ) ae 
a eae 


i=) 


(2.2) 


where y?,, 18 the 100(1 — «)th percentile of the y? distribution. Corresponding equations 
for truncated testing are simular. 

If the r, are unknown, the following methods estimate the values for r, from the data. 
One method uses the likelihood ratio estimate for 7, The second method uses a jacknifed 
version of the first method. The two resulting confidence limits Res and R,, 4.) wie 
and without jacknifing r, respectively. are compared for relative accuracy. Quenouille 
[Ref. 15: pp. 353-360] first reported a method for estimating ratios that reduced the bias 
without increasing the variance. Miller [Ref. 16: pp. 1-15] gives an excellent réviewsramune 
jacknife method which includes a discussion on the application of jacknifing to estimat- 


ing ratios. 


A. LOWER CONFIDENCE LIMIT R,, WITHOUT JACKNIFING 


; , are Cae 
In this case, the maximum likelihood estimate of the ratio v, = 7a 1S 


“Mm 





‘5 fj 
Soe A (28) 
/m 
Where 4,=4H/7, and 4, = max (/,,e¢+¢,4,). The resulting approximate upper confi- 
dence linut for /, is 
y? 
& 1. 2n 
| a 
A m,U(a) i (2.4) 
ee 
ae 
i=] 
k 
Where. 9/0 = erareciic 
=) 
n, 
T,= 7) fy § = V2 oee,k 2a 
j=l 


and 7,, denotes the failure time of the jth test for component i. 


The resulting approximate confidence bound R,,,,) is given by 


( . ke . 1 
Ks tay = aD) — Amcor) a (2.6) 
=I 


B. LOWER CONFIDENCE LIMIT Ree WITH JACKNIFING 
The definitions for ye ,r,,”,, and 7, in Section A are also used in this section. Let 


Aq) denote the estimate for 7, by removing 7, from the data ; iie., 


Pe ica! , | 
ii an OO sie oe) (27) 
Ti 
i=] 
and 
A Hy, a ] ; 4 : ) 
4 mi(o) He ey, — Pe Pes KP LIA Ny, lx] (233) 
ier 


+ eA VES) 
i> hie = : 
i 
where Hw = min(n,, 1,,°¢¢ , 7,} and 
A 
Ass 
A ie) : 
oa ») 
Vile) ‘A J = Le Ly eee iH (2.10) 
* mils) 
ByOredeiine /:,, bY 
y- 
Ae CSP ia 
a ne a i ae 5 
is U(2) — ; (2) 


A* 
2 rT; 
ad 


i=] 


ine corresponding confidence bound R-.,,,, 1s given by 
S S,L(3) = > 


ra” 


k 


A“ Nx * 
Rs ivy = ia =n ryt (232) 
j=] 

When the 7, ¢ = 1,2,+¢++,4 differ considerably, this jacknife estimation proce- 
dure can be unbalanced. That 1s, the number of data points used to compute the jacknife 
estimate will differ from one component to another. It was decided to use this jacknifing 
procedure and determine the effect of differing sample sizes by examining the results of 
the simulations. 

In equation (2.9), the jacknife estimate 1s constructed by using only the first ”° ob- 


A 


servations from each component to obtain ys, 


yey Where mw = min{n,. °°, Hoe 


course it is rather arbitrary to take the first m’ observations from x, , since any ’ of the 


collection of #, values could be used. 


C. SIMULATIONS AND RESULTS FOR EXPONENTIAL CASE 
1. Simulation 
a. Simulation language and package 
The programming language used to simulate this problem was VS 
FORTRAN on an IBM 3033. Also LEXPN in LLRANII was used to generate observed 
exponential random variates. SHSORT was used to perform the sort routines. 
b. Cases and Input parameters 
The six input parameters below determine the conditions for exercising the 
simulation runs. System relability is determined by five parameters. Values of these 
parameters are given in the tables that show the results. The test plan simulated was to 
lest 7 ates UMelealiran, 


* umber Of Component tieesen seeders 


e svstem rehability, A, HOLS) choral 10) 5 )7/s 

¢ significance level, o 0:05 and 

* component time, J, ; small to large (see tables) 
¢ sample size, 1, : small to large (see tables) 


k 

Reliability of a series system 1s expressed as Rs = exp{— >) A,t,}.. We chose arbitrarily 
il 

to determine the failure rates, 4, , from this equation by assuming all /,, to be equal, 


consequent 


c. Replications 
The procedure was replicated 1000 times for each case to get 1000 values 
of Ry,.) and 1000 values of Ree . We order each set of the 1000 values of R,,,,. and 


A 


fein ascends order. Then the two approximate confidence bounds. Kg,,,) y000-.) and 
ae ooo-sy » fOr Ry are the 1000(1 — a)th order values of these two sets of data. Finally, 
for each of the two ordered sets of data, we find the order indices j, and j, for which 
Rs rayj,) and Ree are closest to R, for their respective sets. Then j, / 1000 and 
J,/ 1000 are called the two corresponding simulated true confidence levels. 
2. Results 

Tables 1 through 3 show the results of the 3 cases simulated. The results ind- 
cate that this interval estrmation method using estimates of failure rate ratios will vield 
quite accurate lower confidence hmits for svstem reliability when components have un- 
known constant failure rates. The jacknife method also yields very accurate confidence 
limits which are slightly conservative. That is the 100(1 — «) percentile points of Ree 
Saye i) the tables, are slightly less than the true value of Rg. Consquently, 


ae, = K5)>1—oa . Alternatively, one can say A,;,, is a conservative 100(1 — a) 


ewent rower Coniidence hit procedure for R, . 


Table 1. RELIABILITY OF A SERIES SYSTEM WITH SMALL NUMBER 
( LESS THAN 10) OF SAMPLE SIZES - EXPONENTIAL CASE 


. Lower Confidence True Confidence 
Sumberol | Re cliability Limit Limit 


Compo- ; x 
Re ane of System WO WITH WO WITH 
Jacknifing | Jacknifing | Jacknifing | Jacknifing 
| 2 | 0.9017 | 0.8904 0.9100 
05 | 0.8985 | 0.8889 0.9780 


3 oars Poste ossa [000 
[0s | ostae [09723 | 0.9e00 | 0.9800 


Sample sizesdiem > components are 74 (Umer 
Sample sizesulorsis CORLISS oe (mee | = oi Slag ete Olay abs Mate 
Component times 7, are 2, 5< 7. 8. [O50 7 ee 





Table 2. RELIABILITY OF A SERIES SYSTEM WITH MEDIUM NUMBER 
( LESS THAN 30 ) OF SAMPLE SIZES - EXPONENTIAL CASE 


tea Lower Confidence True Confidence 
Sumber Ol | Reliability Limit ma 


Compo- 
ee of System WEE WITII WO WITH 
Jacknifing | Jacknifing | Jacknifing | Jacknifing 
| 2 | 0.8999 | 0.8978 | 0.8020 0.8660 
[05 [as090 [08083 


hi a ST 
Tos ose | osss1 | asa0 [ost 
209751 [07380790] 0.9000 
Tos [09730 [ora | 09500 [0.9800 


Saniple sizes for S components ane cee loo 
Saple sizes for 15 components ane 347.010, 15, 20815. (0s 20 0a le eee 





Table 3. RELIABILITY OF A SERIES SYSTEM WITH LARGE NUMBER 
( LESS THAN 100 ) OF SAMPLE SIZES - EXPONENTIAL CASE 


ane f Lower Confidence True Confidence 
me! Reliability Limit Limit 


ils 
achuvae came Jacknifing | Jacknifing 
eg oe ae 
| 8 acpi Caos Poe 

0.9520 
Se 
A 


Sample sizes for 5 components are 30, 63. 75, 98. 26 
eemmeeie sizes for 15 coniponents are 15, 40, 33. 17, 26, 67, 50, 65, $0, 32, 95, 100, 15, 
45, 30 





Il, INPFERVAL ESTIMATION PROCEDURE - WEIBULL CASE 


A. LOWER CONFIDENCE LIMIT R,, 
Consider a series svstem with A coniponents. Let the time to failure X, of comme: 


nent / have a Heibull distribution with density 


ae: = : 
fil) = AE BaP" exp Gin}. 4 >. (3.1) 
Then 
R(t) = exp — (4,0)"'4 > 0 (3.2) 
Det AP hs y i ; 
and 
y i= 
Rany= exo) = >) Hit | me ~ind) ntl Oo exe (Jee) 
i=] el 
where A =/% ,/, = max/ and 7, = ///.. Ifthe 6, are known, A” has cote 


failure rate 7% and the procedures 11 Chapter JI can be used to obtain /(2yeemne 
Reig) With [, replaced by 7? in equation (2.5). itis well kKitiewm that 7} = lives aeaes 


i 


an extreme value distribution with CDF 


as 
Fy) = I- exp) — ora 
l 5 
where ¢ = inl? jwand 6 —= 3) 6. 
Engelhardt and Bain (Ref. 17: p. 323] have developed the following simple unbiased 


estimators for €, and 6, . using ordered values 1}, =1nX,, 


S n 
: 5 ; 
= yy n—S » Vy 
] j=i j=S+) oa 
B. = nik Om 
et 


ae 


where s = (0.841,] = largest integer < 0.84 n, and X,,. 1s the jth order statistics from the 


SQMIpIC Gi Size as, Oly ane SISO 


SX > 
I 
—— 
2 
a 
iD |— 
ee” 
| 
So 
+ 
Gr 
— 
oe) 
an 
eee” 


miiere y=0.57/2 and +,=S Y,/n,. Let 


- 


J=1 


and 


1 a aia) |e) 


i 


i aes kK. (oa) 


We approximate the distribution of 7, by the exponential distribution with failure rate 


28 = A> and proceed as in Chapter I]. Define 


i 


+ Nn; (3 S) 
eee a 
; 
ni A A 
mane 7, = +f, § = 1.2,¢++,k. Let 4,=max/, and 
=} { 
Ne 
A Zj 4 
n> one (3.9) 
ip 


7 Xx 2n i 
2 a nT; 


and the corresponding approximate lower confidence limit As. for A(z) 1s given by 


k 
Ks 143) = exp | ae Ame), rab (321 1) 


| 


! 


k 
ance 7 = > 7... 
a 
This procedure is labeled the Formula procedure in the tables that follow in this 
section. Its distinguishing feature is the equation for f, given by equation (3.6). An al- 


ternative procedure, labeled the Newton - Raphson procedure, estimates f, using the 


maxunum likelihood procedure which 1s obtained using « Newton - Raphson approxi- 
mation method. Equations for 8, under the Newton - Raphson procedure are provided 


in Appendix A. 


B. SIMULATIONS AND RESULTS 
1. Simulation 
a. Simulation language and package 
The programming language used to simulate this problem was VS 
FORTRAN on an IBM 3033. Also LEXPN in LLRANII was used to generate observed 
exponential random vaniates. SHSORT was used to perform the sort routines. 
b. Cases and Input parameters 
The six input parameters below determine the conditions for exercising the 
simulation runs. System rehabilitv is determined by five parameters. Values of these 
parameters are given in the tables that show the results. The test plan simulated was to 
test-72sikelns Witla taal. 


¢ number of component tvpes, A: 5 and 15 


e system reliability, R, NOE! sesh WES TS) 
e significance level, « OS and. Ue 
¢ component time, /, ; small to large (see tables) 
¢ sample size, 7, : small to large (see tables) 
k 
Reliability of a series system is expressed as R, = exp{— S42}. Failure rates /, can 


t= 


be determined from that equation by assuming all 4,7, to be equal, consequently, 


— Ihara. 
7 =e 
c. Replications 

The procedure was replicated 1000 times for each case to get 1000 values 
of Kop. calid COU eames aii Rete) . We order each set of the 1000 values of K,, ama 
a in ascending order. Then the two approximate confidence bounds, K. 7.) joccqesyeae 
Rev akitane: , for R, are the 1000(1 — o)th order values of these two sets of data. Finally, 
for each of the two ordered sets of data, we find the order indices j, and J, for which 


Rsreyjy) ANG Re roy;,) are the closest to Rs for their respective sets. Then j, / 1000 and 


J, / 1000 are called the two corresponding true confidence levels. 


2. Results 
Tables 4 and 5 displav the results of the simulations and determine the accuracy of 


Regie) given in equation (3.11) as a lower confidence limit procedure for svstem rehability, 


R,, for parameter values /,, 8, 4, and R, given in the tables. The terms in the tables have 


t 


na 


the same meaning as the corresponding terms in Tables | through 3 which were ex- 
plained in Section 2.C.2 . The procedure would be exact for the Formula method if the 
values in the Formula column equal the corresponding numbers in the same row in the 
Reliability of Svstem column. In Table 4 for example, the 80th percentile point of .9205 
for the Newton - Raphson procedure is more accurate than the Formula procedure 
which has an 80th percentile point of .9324. The last column shows, however, that what 
we have called an 80% lower confidence limit procedure is in fact closer to a 61% pro- 
cedure. The accuracy improves if the sample size 1s increased from 15 to 30 as indicated 
in Table 5. The accuracy improves even more if the number of components in the svstem 
increases from 5 to 15 as indicated in Tables 4 and 5. 

The results in Tables 4 and 5 show that the Rie method given by equation 
(3.11) is too optimistic - especially for small sample sizes and for systems with a small 
number of components. Modifications to this procedure are needed. These tables also 


indicate that the Newton - Raphson method 1s superior to the Formula method. 


Table 4. ACCURACY OF R,, AS A LOWER CONFIDENCE LIMIT FOR R; 
WITH SAMPLE SIZES OF 15 - WEIBU! L CASE 
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Table 5. ACCURACY OF R,, AS A LOWER CONFIDENCE LIMIT FOR R, 
WITH SAMPLE SIZES OF 30 - WEIBULL CASE 


x Lower ie True Confidence 
Number of Reliability Limit Limit 
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P p 
| 


Initial beta for 5 components are 2.2, Fe 2:6G2-00 225 
Initial beta for 1S components are 2.2, 7-4, 2.6, 7302 532 3 2, 2.2 20 ee 
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IV. INTERVAL ESTIMATION PROCEDURE - NORMAL CASE 


A. BACKGROUND 
In recent years it has become popular to model mechanical rehabilitv as PLY > Y) where 
X denotes strength and Y denotes stress. This formulation of reliability is an average of 


P(Y>y| ¥=3), since 
P(X > Y)= Ey P(X >y)= | PXY>Y)LO) ay. 


An alternative model is the worst case approach. In the worst case model, one se- 
lects a worst case value of }, savy, . and then designs the strength. Y ,of the component 
so that PLY > 4,) = R, where R, is a rehability requirement. If VY has a normal distrib- 
ution then this requirement imposes constraints on the mean and variance of Y . This 


1s usually done in a manner to comply with standard “safety factor” procedures. 


The average model, P(X > }), uses two random variables and 1s subject to more 
Fandom error than the worst case model. The accuracy of the expression P(X > }) for 
values of this expression close to one 1s highly suspect due to deviations in the tail 
probabihties of both X and ¥ from those assumed in the model. It is common to assume 
that both Y and Y have normal distributions. Truncated normal distributions would be 
more appropriate for many types of mechanical equipment. Harris and Soms [Ref. 18: 
pp. 650-663] discuss implications of this problem. Very significant errors in point and 
interval estimation for reliabilitv are readily demonstrated when VY is truncated normal 
but assumed to be normal in the more simple model which specifies that R = P(X > ¥). 
Table 10 shows this effect when YX is truncated above at u + 1.6450 , where Z, is the 100 
(l—«)th percentile point of the standard normal distribution. Church and Harris [Ref. 
12: pp. 49-54] and Downton [Ref. 13: pp. 551-558] have developed approximate confi- 
dence intervals for PLY > ¥) when X 1s normal with unknown mean and variance and Y 


has the standard normal distribution. 


B. EXACT AND APPROXIMATE INTERVAL ESTIMATES 
Minimum variance unbiased estimators (MVUE) for R= P(Y > y) are well known 
when X 1s normal with unknown mean and known variance and also when the variance 


is unknown. In the former case, Lieberman and Resnikoff [Ref. 8] developed the result 


A ve) 
kno ee) 
ox / n— | 


H 


Which is MVUE for R where @ 1s the standard normal cumulative distribution function. 
When the variance is unknown, several versions of an integral expression for P(\ > y) 
have been developed by Lieberman and Resnikoff [Ref. 8] , Basu [Ref. 19] and Folks and 
others [Ref. 9]. 

Exact lower confidence interval estimates for PLY > +) when X is normal with un- 
known mean and variance involve the non-central ¢ - distribution. Owen and Ilua [Ref. 
11] developed tables of the lower 90% and 95% confidence limit values R, for PUX ay) 
based on the non-central ¢ - distribution. These values are tabled for values of & in the 
range -3.0(.2)6.0 and sample sizes n = 2 (1) 18, “21 (3) 30. 40 (20) 100 Fxvieie 
A = (x —3)/s and x and s are the sample mean and standard deviation. Their tables are 
essentially exact. An approximation to their exact tabled values is given bv R; where 


K, = ages 


: 2 0 


ine Te 


(4.1) 


4,n—-| 


kK=(xX—y)/s andi 


grees of freedom. Equation (4.1) was developed in this thesis. It 1s an extensive ad hoc 


.n-) 2S the 100(] — «)th percentile of the 7 distribution with » — 1 de- 
modification of an equation developed by Church and Harris [Ref. 12]. Tables 6 and 7 
display values of Rk, = @M(y),y,7', R, . and KR, —R, for &k = |. 2, 3, 4, sample sizes x 
= 10, 18, 30 and confidence levels 90% and 95% . K, and y denote the exact owas 
confidence linuts and corresponding ®7!(R,) values from Owen and Hua [Ref. 11]. Both 
y’ and R; are given by equation (4.1). The accuracy of the approximate confidence in- 


terval 1s quite ood relative tol tite values tor A. oem by Omen adel asi comme 


16 


Table 6. APPROXIMATE ( Ry ) 
LENRS FOR P(X SF) 


13 | 1.3334 


AND EXACT ( R, ) 90% CONFIDENCE 
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Table 7. APPROXIMATE ( RAR? ) AND EXACT ( R, ) 95% CONFIDENCE 
Evils FOR PY = 17) 


Crees: Posssi | eam [0.60913 | 003140 

3 | 0.95565 0.93685 | -0,01910 

2.07743 0.98075 -0.00956 

Tr Focasns | o.stonT | 0.9202 | 0.6ss96 [0.00605 
1.30221 0.89777 | -0.00571 

IS | 1.7396 oe : ve 
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Tables $ and 9 display the results of computer simulations with 1,000 replications 
to check the accuracy of the R; method for 80% and 90% lower confidence limits for 
PNY>y)=R fory = 3. with various values of o and uw determined so that R equals 
the values shown. The procedure would be exact at the 80% level if the values in the 
column labelled o = .2 equal the corresponding values of R in the same row. The “true” 
confidence level corresponds to the index i(R) of 1,000 ordered values of R;,,, for which 
Riear = & - For example, in the seventh now of lables; Ro — .950), co, = 20. Neen 
Riva s0 WAS -9575, and Ky cipyoo WAS -9533. ASO Rico -cg = -950 and Kj igen, eee 
Tables § and 9 indicate that the R; procedure given by equation (4.1) 1s quite accurate 


at the 80°% level of confidence for the cases simulated. 
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Table 8. ACCURACY ANALYSIS OF Rx PROCEDURE FOR 80% AND 90% 
LOWER CONFIDENCE LIMITS FOR P(\X>3) WHEN X_ IS 
NORMALLY DISTRIBUTED 
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Table 9. ACCURACY ANALYSIS OF RF PROCEDURE FOR 80% AND 90% 
LOWER CONFIDENCE LIMITS FOR PLY >30) WHEN X_ IS 
NORMALLY DISTRIBUTED 
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Table 10 displays the inaccuracies of the R; lower confidence limit procedure for 
P(X > 3) using equation (4.1), which assumes X is .\(j. o7), when in fact X has the dis- 
tribution of a normal random variable with mean p and variance o? that has been trun- 
Gied at + |(643¢. Note that p+ 1.6856 > > 3, becalise P(X > 3)=.95,.99 and .995. 


Examination of Table 10 reveals gross inaccuracies; consequently, even when the 
distribution of Vis truncated far above the value y, the exact lower confidence limit for 
P(X > y) can be greatly in error when computed under the assumption that \V is normal. 
This problem will be compounded when one is computing “exact” confidence intervals 
for P(Y > ¥) assuming both V and FY are normal when in fact one or both may be trun- 
cated normal. This suggests that modeling mechanical rehabilitv as PLY > }) mav be 
more risky than more standard engineering approachs for modeling mechanical reliabil- 


itv which include the notion of safety factors. 


Table 10. ACCURACY ANALYSIS OF R7 PROCEDURE FOR 80% AND 90% 
LOWER CONFIDENCE LIMITS FOR PX > 3) WHEN X IS TRUN- 
CATED NORMAL 
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V. CONCLUSIONS AND RECOMMENDATIONS 


The lower interval estimation procedure for reliability of coherent systems which 
was developed in Chapter If appears to be accurate, easy to use, and applicable to co- 
herent systems. Although this procedure assumes that failure times of all components 
of the svstem are independent and have exponential probability distributions, it can be 
easily extended to systems with component failure distributions that can be transformed 
into exponential failure distributions ; e.g.. HWerbull distribution with known shape pa- 
rameter. This procedure has potential for beng combined with a similar procedure for 
svstems with cyclical components. The combined procedure would provide for use of 
binomial component test data and exponential component test data to compute lower 
confidence limits on the reliability of coherent systems with both cyclic and continuous 
components. Such an extension could use a failure rate ratio estimation procedure sim- 
ilar to that developed m this thests. Such a method should be explored. 

The interval estimation method for the rehabilitv of a system with components that 
have Hefbull failure trmes 1s not sufficiently accurate to be applied to systems that have 
10 or fewer components each with ten or fewer tests. This procedure needs further study 
and refinement. 

The approximate lower confidence linut for component relrability PLY > 1) when 
is normally distributed with unknown mean and variance 1s very accurate. It does not 
require an extensive set of tables such as those developed by Owen and Hua[Ref. 11]. 


but only requires the use of the standard 7 tables. 
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APPENDIX A. MLE OF WEIBULL PARAMETERS BY THE NEWTON - 
RAPHSON METHOD 


Let X ~ HWEl(a, 6). Then the likelihood function for the first r ordered observations 


from a random sample of size 7 1s given by 


Fi * 0% «rad = a Ts a fo = Fi 
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ELEY oof [Yeo] 


i=} 





i= 


Setting the partial derivatives of this likelihood with respect to « and f equal to zero 


gives the MLE’s & and f as solutions to the equations 





A = 
a i 
i=} 
and 
nN 
al 
», xP “A, 
a H 
Where #7 = r. It can be shown that these equations have unique solutions which are the 


maximum likelihood estimates. The NEWTON - RAPHSON method for solving an 
equation e(f)=0 1s to determine successive approximations f, where 


Ne e(B.)/g'(B). Therefore, the estimates of « and £ can be solved by letting 
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APPENDIX B. FORTRAN CODE FOR INTERVAL ESTIMATION 
PROCEDURE - EXPONENTIAL CASE 
PROGRAM EXPONE 


Sevcvevevesodesledededvdevededeslededescdkcaedededcsdcckcsevledesededcsevedesleslestesicnevesedesicttesledeslcseslesicvesledesicseslessscslevedcdesvededte 


se * 
¥ THIS IS A PROGRAM TO COMPUTE THE TRUE CONFIDENCE LIMIT OF * 
% A SERIES SYSTEM WITH / WITHOUT JACKNIFE METHOD AND * 
ve COMPARE THE DIFFERENCE OF THOSE KESUL TS: ¥ 


Ce res Cred cre) ce) CC) Cnr) Ce Ce) Cn CC C) Cr C 2) C2) C).C) OC) Gl OC) Cr) 9 I I a oa a a ae 


se 


ae 


N 

RSYS 
ALPHA 
Ah IUialls 
SAMPLE 


B 
BIGLAM 
CASE 
CHISQR 
CLR1 
CLR2 
CTIME 
Be 
DEE 
LAMBDA 
LAMHAT 
LAMHM 
LAMHST 
LAMMAX 
LOTT 


ROMIT 
ROMSUM 
RSTAR 
RVAL1 
RVAL2 
Ah 
ZALPHA 


we we Wwe we Wwe we we wee Wt we We we we et ee ee we we 


we we Wwe Wwe we We Wwe Wwe We He We 


BELOWS ARE GIVEN OR ASSUMED. 


NUMBEK “OF (COMPONENTS 2 UyYer 

RELIABILITY OF A SERIES SYSTEM 
SIGNIFICANCE LEVEL 

TEST TIMe 

SANPLE SIZE FCkK EACH ty EhawheCONPONERT 


THESE ARE VARIABLES USED. 


TEMPORARY ARRAY FOR EXPONENTIAL RANDOM VARIATE ?* 


LARGEST VALUE OF LAMBDA 
COUNTER 
CHI-SQUARE VALUE FOR GIVEN ALPHA AND DF 


TRUE CONFIDENCE LIMIT FROM NON-JACKNIFING 


TRUE CONFIDENCE LIMIT FROM JACKNIFING 
TEST TIME OF INDIVIDUAL COMPONENT 
DEGKEE VOR PE REEDOMSHORS Chl SQUARE: 
ABSOLUTE, VALUE“ OR DYGReRGNGr eer TWEEN kl 
FAILURE RATE OF EACH COMPONENT TYPE 
LAMBDA HAT FOR JACKNIFING 

FINAL LAMBDA HAT FROM NON-JACKNIFING 
FINAL LAMBDA HAT FROM JACKNIFING 
LARGEST VALUE OF LAMHAT 

LAMBDA WITH 1 COMPONENT OMITTED 
TEMPORARY ARRAY 

TEMPORARY ARRAY 

INITIALLY COMPUTED RELIABILITY FOR EACH 
COMPONENT TYPE 

COMPUTED RELIABILITY BY NON-JACKNIFING 
COMPUTED KELIABIMI iy SEY SJACENIEING 
RATIO OF LAMBDA FOR NON-JACKNIFING 
RATIO OF LAMBDA FOR JACKNIFING 

RATIO WITH 1 COMPONENT OMITTED 

USED FOR GACKNIFE. SUMO aah Ori i 1 
FINALLY COMPUTED R VALUE ( R STAR ) 

RC 500 * (1-ALPHA) ) FOR NON-JACKNIFING 
R¢ 500 * €1-ALPHA) ) FOR JACKNIFING 
TOTAL TEST Tih OF EACH COLO rN ie er 


& KZ 


RIGHT PERCENTILE POINTCNORMAL DISTRIBUTION) 


Cl Cl (acca 


10 


20 


25 


sc 


Jesesedececicccccdevevevesesestcstcacsccdededededededescacdcvdedevestslescvedesededeskdedcdcdedededesedevevededededededededededevevedvede ve 


PARAMETER (NN = 500) 


REAL TIME(15), ALPHA(2), ZALPHA(2), RSYS(2), LAMBDA(15), P 
REAL LAMMAX, LAMHAT(15), B(100), RSTAR(15), LAMHM 

REAL RATIO1(15), RATIO2(15), BIGLAM, DIFF(15) 

REAL LOMIT(15,100), CTIME(15,100), ROMIT(15,100), T(15) 
REAL ROMSUM , SUM1, SUM2, LAMHST, R1(NN), R2(NN), KEY1(NN) 
REAL RVAL1(15), RVAL2(15), CLR1(15), CLR2(15), KEY2(NN) 
REAL DIFFR1(15), DIFFR2(15), R(15) 


INTEGER SAMPLE(15), N(2), DF, CASE 


DATA N / > 0 / 

DATA RSYS /. >» mO7 5 of 

DATA TIME ei 8 ONee 1O wie eG aie a 

DATA SAMPLE / 30,63,75,98,26,15,7,5,30,20,10,7,13,21,30 / 
DATA an Oo Osta 

ims 1SEED) / 1736 / 


/** SUBROUTINE ZTABLE COMPUTES THE RIGHT PERCENT POINT ZALPHA 
/** FROM RIGHT CUMULATIVE PROBABILITY ALPHA 
POmlOnl = 1.02 
CALL ZTABLE(ALPHA(I), ZALPHA(I)) 
CONTINUE 
CASE = 0 
/* II IS INDEX FOR N 
DO 150 II = 1, 2 


/* COMPUTE THE DEGREE OF FREEDOM FOR CHI-SQUARE 


DF = 0 
DO 20 I = 1, N(II) 

DF = DF + SAMPLE(1) 
CONTINUE 
De 2 DP 


/* FINDING COMPONENT TYPE THAT HAS THE MINIMUM NUMBER 
/ -SOPVSANPEE SIZE 


NSTAR = 999 
DOR Sale le NGL) 
TE CSAMPLECT) 2 bE. INSTAK) THEN 
NSTAR = SAMPLE(1) 
ENDIF 
CONTINUE 


% / 


/* JJ IS INDEX FOR RSYS se 
DOL IaG — Taee 


/* COMPUTE LAMBDA FROM THE GIVEN EQUATION AND ¥ 
/* FIND THE BIGGEST LAMBDA a / 


BIGLAM = 0 

DO 30 K = 1, NCII) 
LAMBDA(K) = ( -ALOGCRSYS(JJ)) / NCII)) / TIMECK) 
IFC LAMBDACK) .GE. BIGLAM) BIGLAM = LAMBDA(K) 
R(K) = EXP( - LAMBDACK) * TIME(K) ) 

CONTINUE 


/* COMPUTE CHI-SQUARE VALUE FOR GIVEN ALPHA AND DF ie) 


/* LL IS INDEX FOR ALPHA * / 
DO 130 LL = 1, 2 


IF ( DE .EO. 1) THEN 
P = ALPHA(LL) / 2 
CALL ZTABLE(P, ZALPHA(LL)) 
CHISQR = ZALPHA(LL) *** 2 


ELSBSUE we 0P EO 2 ies 
CHISQR = -2 * ALOG(CALPHACLL) ) 


ELS Ear GbE or: so ene 
GCHISOR +=) 8) "> Cll= 24 (os: DE ae 


ZALPHAC( LL) * SQRT( 2./(9 * DF) ) )ijagaee 
END 


CASE = CASE + 1 
DIFF(CASE) = 0 


DO 120 L = 1, NN 
LAMMAX = -99. 
fe Sa ak, NCS 
/* GENERATE EXPONENTIAL RANDOM ve / 
/* NUMBERS WITH MU = 1 ‘ 


CALI LEXPNCTSEED er cain ( Lyi, 0) 
TCI) = 0 


DO 40 J = 1, SAMPLE(T) 


/** CONVERT TO EXPONENTIAL RANDOM */ 
/** NUMBERS WITH MU = LAMBDA AND */ 


C 
C 
40 
C 
C 
50 
C 
C 
60 
C 
C 
70 
80 
70 


/* ADD THOSE FOR EACH COMPONENT 


/* TYPE 


B(J) = B(J) / LAMBDA(I) 

T(I) = T(I) + BCJ) 

CTINEC Io) = BC J) 
CONTINUE 


LAMHAT( I) = SAMPLECI) / TCI) 


/* FINDING MAXIMUM LAMBDA HAT AND ITS 


/* INDEX 

IF ( LAMHAT(I) .GE. LAMMAX ) THEN 
M=I 
LAMMAX = LAMHAT(I) 

ENDIF 


CONTINUE 


/* RATIO1 IS FOR WITHOUT JACKNIFE 
/* RATIOZ TS FOR Wii JACKNIFE 


DO 60 I = 1, NCII) 


RATIOI( 1) = LAMBDACI) / BIGLAM 
RATIO2(1) = LAMHAT(I) / LAMMAX 
CONTINUE 


/* PART OF JACKNIFE METHOD FOR LAMBDA 
/* WITH OMIT 1 VARIABLE EACH TIME 


DOeSO T= ee el. ) 
DOr 307 = SAMPLECT) 


LOMtl (iJ) = 0. 
DO 70 K = 1, SAMPLE(T) 


IF (J .NE. K) THEN 


7 / 


ay 7 
vc i 


DY / 
a / 


re 
“ey pay 
— 


Or Mee. tI) = LOMITCI,J)-+ CTiIMEC LK) 


ENDIF 


CONTINUE 


LOMIT(I,J) = (SAMPLE(I)-1) / LOMIT(I,J) 


CONTINUE 
CONTINUE 


Zo 


i i ie el 


C2G7- Cac] 


/* ADAPT ABOVE RESULT TO OUR EQUATION TO */ 
/* GET THE RELIABILITY AND TRUE CONFIDENCE */ 
/* LIMIT */ 
SUM1 = 0. 
SUM2 = 0. 
SUM3 = 0. 
SUM4 = 0. 
DO 110 I = i, NCIT) 
/* NON JACKNIFING ( ORIGINAL ) * / 
SUM1 = SUM1 + RATIO1(I) * T(I) 
SUM2 = SUM2 + RATIOI(I) * TIME(I) 
/* WITH JACKNIFING ve / 
ROMSUM = 0. 
DO 100 J = 1, NSTAR 
ROMIT(I,J) = LOMIT(I,J) / LOMIT(M,J) 
ROMSUM = ROMSUM + ROMIT(I,J) 
100 CONTINUE 
RSTAR(I) = NSTAR * RATIO2(I) - 
3 (NSTAR - 1) * ROMSUM / NSTAR 
SUM3 = SUM3 + RSTAR(I) * T(I) 
SUM4 = SUM4 + RSTAR(I) * TIME(I) 
110 CONTINUE 
/* R1 ; RELIABILITY OF A SYSTEM WITHOUT * / 
Von JACKNIFE ze / 
/** R2 ; RELIABILITY OF A SYSTEM WITH te / 
/* JACKNIFE / 
LAMHM = CHISQR / (2 * SUM1) 
R1(L) = EXP( - LAMHM * SUM2 ) 
LAMHST = CHISQR / (2 * SUM3) 
R2(L) = EXP( - LAMHST * SUM4 ) 
IF ( ABS(R1(L) - R2(L)) .GE. DIFF(CASE) ) THEN 
DIFF(CASE) = ABS(R1(L) - R2(L)) 
DIFFR1(CASE) = R1(L) 
DIFFR2(CASE) = R2(L) 
ENDIF 
120 CONTINUE 


/** NONIMSL LIBRARY 'SHSORT' WILL SORT R1, R2 
By SHELL SORT ALGORITHM 


/* 
AY 


ae / 
ate 
PAY 


30 


CALL SHSORT(R1, KEY1, NN) 
CALL SHSORTCRZ, KEY2Z5°sN) 


MM = NN * (1 - ALPHACLL)) 

RVALICCASE) = R1(MM) 

RVAL2(CASE) = R2(MM) 

/* SUBROUTINES IND FINDS hia wINDES OF Ri. K2 
/* WHIGH THE VALUE OF 171 IS CLOSEST TO RSyYS. 


CALL FINDJ(R1, NN, RSYS(JJ), J1) 
CALL FINDJ(R2, NN, RSYS(JJ), J2) 


CLR1(CASE) = J1 / FLOATCNN) 
CLR2(CASE) = J2 / FLOAT(NN) 
130 CONTINUE 
140 CONTINUE 
150 CONTINUE 
WRITE( 6,600) 
WRITE(6,650) 
WRITE(6,670) 
CASE = 1 
DooetO Il-= 41, 2 
DONO 0 = 1.2 
Doo On i= ad 
WRITE(6,700) N(II), RSYS(JJ), ALPHA(LL), RVAL1(CASE), 
ve RVAL2(CASE), DIFF(CASE), CLR1(CASE), CLR2(CASE) 
WRITE(6,888) DIFFR1(CASE), DIFFR2(CASE) 
CASE = CASE + 1 
190 CONTINUE 


200 CONTINUE 
WRIgeiCo.7 77) © SAMPLECS).-J=1,NCEI) ) 
Whee C6. 992) (°kUJ), J=l,NCl1) ) 

210 CONTINUE 


600 FORMAT('1' ,5(/),7X, '**** RELIABILITY OF SERIES SYSTEM *?7' ) 
650 FORMATC///;, T50, Rl, CL1 ; WITHOUT JACKNIFING’ , 
PSO. Ro eC? ; WITH JACKNIFING' ) 


670 FORMAT(///, T6, NUMBER OF' ,T19, RELIABILITY’ ,T33, ‘ALPHA’ ,T46,'R1' 


T56,' R2' Mies. “MAX (RI = RZ). ,T79, TRUE’ ,T87, TRUE ,/, 
T6, COMPONENTS ' elon @ SYSTEM’ ,T79,'C. L. dolore Ga. 


/,T5,89('-")) 


9 ' 


ad 
% / 


bd 


700 FORMAT( /,T8,15,T22,F5. 3,T33,F5. 2,T43,F7.4,T53,F7.4,T65,F8.4,1778, 
xc F7.4,T86,F7.4) 


777. FORMAT(/,T3,'SAMPLE SIZE ARE ', 5(2X,13) ) 
888 FORMAT(T60,'Rl=',F6.4,' R2=' ,F6é.4) 
999 FORMAT(/,T3, RC1) 3) .o0 bese) 


STOP 
END 
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SUBROUTINE ZTABLE( ALPHA, ZALPHA) 


<< SUBROUTINE ZTABL2 COMPUTES RIGHT PERCENT POINT ZALPHA PROMS 
<< RIGHT CUMULATIVE PROBABILITY ALPHA >> 


C2 2G) C2 


REAL ALPHA, ZALPHA 
IF CCALPHA ~GT. 0.0) .ORS)( ALPRAS2 ER ses pe ien: 


W = - ALOG(4 * ALPHA * ( 1 - ALPHA )) 
ZALPHA = SQRT(W * (2.06118 - ( 5.72622 / (W + 11.6406)))) 


IF (ALPHA . GOT: 2025) AEN 
ZALPHA = - ZALPHA 
NBN S 
ENDIF 


BREGURN 
END 


ale alantsateulontantaaisetaatsats alton euloute eu! Seseslevesledclevevevededeveve te sioclontontc fa wlowtewtoatorutsuleatsulerteutoutoutwulactoutsutsatontoutsulwatsuteoule Jerlseleses'es'eve wleute ale wleateuts 
OR Oe CR ex ee a FX FY CL EY 4 EX EN raver is river ty VERNER REL 4d EX FL CX UH Ed EC rier iy vee 8 eu as ee rier ty ey @e eh eR 4% raver ts ever rkeris @% TITAS @e eh 4B er eh @% SR EN FR FD ER ER ER EN CR OY Cd FR 
ae ae 
é @. 


ZO ON TREN TN EN FS EN EN ON EN EN OVER ON ER ER EN EN EN GR OO EN EN ER EN SEEDER EC EN ENED GR EN ES ELEN ER EN EU EVER ANON EN EN ON EN EO EN EN ER EC ER ELEN EL EN EN EU EN EW ENED EV ELEN ON ENON OD 


SUBROUTINE SHSORT(A,KEY,N) 
DIMENSION ACN) ,KEY(N) 
M1=1 
6 M1=M1*2 
LF AC Mile e LE ae) GOO hae 
NI=tay Zeal 
MM=MAXO(M1/2,1) 
GOe1O Za 
20 MM=MM/2 
TE CMM LES 0) AGO moar o 
21 K=N-MM 


2 


HU 


100 


DO 1 J=1,K 
LI=J 
IM=II+MM 


TE RGACIMaeGE. AC LEDs GORMer 1 
TEMP=AC IT) 

P= Ki yer ) 

AC II)=AC IM) 
KEY( II )=KEY( IM) 

AC IM)=TEMP 

her CIN j=f1 

= 1 1 -MM 

ill. Cleese GO-TO 11 
CONTINUE 

GO TO 20 

RETURN 

END 


Ievedededslsdedvdevesededekideclvckkkakkkkevedcdledkedcakdedescak ovakocdksleckdeclescakeseckscsicskclcclcstccvcoksdleskckscotvestvesedesesk assess scctese 


ule 
® 


D 
Phy 


Veleveddededdeveddeddesiekkeckeickicdikdvdedckdedcdedekdedeskdevedcvedesoveteledckevevedevedkevcvevedesevenkevevevededevevevedeveves's 


a > Us Ee | 
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SUBROUTINE FINDJ(A, NN, R, J) 


<< SUBROUTINE FINDJ FINDS THE INDEX OF ARRAY A WHICH THE eS 
<< VALUE OF IT IS CLOSEST TO R. >> 


REAL ACNN), R, VALUE 
INTEGER J 


VALUE = ABS(A(NN) - R) 
DO 100 I = NN-1, 1, -1 


TP e@eArSCAC I) = R)o. LI. VALUE) THEN 
VALUE = ABS(ACI) - R) 


ELSE 
ges ler 
RETURN 
ENDL 
CONTINUE 
RETURN 
END 


Ss 


APPENDIX C. FORTRAN CODE FOR INTERVAL ESTIMATION 
PROCEDURE - WEIBULL CASE 
PROGRAM WEIBUL 


Sevesesesesevedevestevedededevededededededededevevededededevededevevtededevedevedeseslcsccestcvesedevedescveseslevlcdesevesesvevestevesedese 


oy THIS PROGRAM COMPUTES THE RELIABILITY OF A SERIES SYSTEM * 


ve USING TWO DIFFERENT METHODS FOR THE CASE WHEN FAILURE ¥ 
ve OF EACH COMPONENTS IS WEIBUGLY Paste iBUrED: % 
“we AND ALSO COMPARES THE ESTIMATE OF SHAPE PARAMETERS ( BEes 
%* WHICH 1S COMPUTED FROM GyVO5Dit EERE) ee anODs, ¥ 
‘c els 
% THE SERIES SYSTEM CONSTRERED INGiMiserOnChAalietsAs ve 
% N COMPONEND-TYPES AND EAGH COMPONENT HAS SANE NUMBERSOE ¥ 
¥ SoberGe olZho. * 


BELOWS ARE GIVEN OR ASSUMED ; fe 


ace Guta cy Gseu Ipc) C2 Cl C)Ciey G).C))C) CC) Glee GC) Cl Ca CIC) ©) CC) 02 C7 OS) 0) 6) Oo a) OO a 


BETA 
ALPHA 
X 
RSYS 


THESE ARE 


BETHAT 
BETNEW 
CHISQR 
CLR 
CLRNEW 
EAPO 
KSUBM 
LAMBDA 
LAMBIG 
LAMHAT 
LAMMAX 
LAMMU 
LAMSTR 
LAMNEW 
LAMUNR 
RATIO 
RATNEW 
RHAT 
RHTNEW 
XBENEW 
XHAT 
XI JHAT 
AIJNEW 
XNEW 
IG O3s UT Ab 
W,WEIB 


e 
bd 
° 
y 
> 
bd 


SHAPE PARAMETER OF WEIBULL DISTRIBUTION 
SIGNIFICANCE LEVEL 

TEST TIME 

RELIABILITY OF A SERTES@s yod Et 


VARIABLES USED ; 


we Wwe te: We Wwe we we we Wwe we wt we we we wee let el Ole et lee lee let Olle lw 


BETA HAT FROM THE FORMULA 

BETA HAT FROM NEWTON - RAPHSON METHOD 
CHI-SQUARE VALUE FOR GIVEN ALPHA AND DF 
TRUE CONPIDENCE, Livi) ERG tine Pore EA 
TRUE CONFIDENCE DIMI PROS ThE N=keMEciop 
EXPONENTIAL RANDOM NUMBERS WITH MU=1 

K SUB M VALUE 

SCALE PARAMETER OF WEIBULL DISTRIBUTION 
MAXIMUM VALUE OF LAMBDA FROM THE N-R METHOD 
ESTIMATES OF LAMBDA FROM THE FORMULA 
MAXIMUM VALUE OF LAMBDA FROM THE FORMULA 
FINALLY COMPUTED LAMBDA FROM THE FORMULA 
LAMBDA TO THE BETA 

ESTIMATES OF LAMBDA FROM THE N-R METHOD 
FINALLY COMPUTED LAMBDA PROM THE N-kK METHOD 
RATIO OF LAMBDA FOR THE FORMULA 

RATIO OF LAMBDA FOR THE N-R METHOD 

RC 5005 G1 -ALPHA) Ok inh rOrtiGs 

R( 9500) * (1 -ALPHA) )ShOh SHE Nak MEdheD 

X TO THE BETNEW 

SUM OF XIJHAT 

W TO THE BETHAT 

W TO THE BETNEW 

SUM OF AIJNEW 

X TO THE BETHAT 

WEIBULL RANDOM NUMBERS 


34 


LP ial GP sp? BOP fe oP es 


OB ie ae Ds a i a 


i ia D1 | 


10 


se 


se 


Seveseses' 


ZALPHA ; RIGHT PERCENTILE POINT (NORMAL DISTRIBUTION) 


ante ate mt 


ssevevevedoveddodcavdedekdevetevesledevesvesesesvevedcvesedesvescsevevesescsesedescdevescvevesessvedesess ved ve sevedese de 


PARAMETER ( N=5, M=30 , NN= 500 ) 
PARAMETER ( N=15, M=30 , NN= 500) 


REAL 
REAL 
REAL 
REAL 
REAL 
REAL 
REAL 
REAL 
REAL 


INTEGER 


DATA 
DATA 
DATA 


DATA 


DATA 


DATA 


DATA 


DATA 
DATA 


DATA 
DATA 
DATA 
DATA 
DATA 
DATA 
DATA 


CASE 


BETA(N), LAMBDA(N), X(N), XHAT(N), XTOBET(N), XX(N,M) 
WCN,M), YCN,M), BETHAT(N), WBHAT(N,M), LAMHAT(N) 
RATIO(N), LAMSTR(N), LAMMAX, RSYS(2), KSUBM(N), LAMMU 
EXPO(M), ALPHA(2), ZALPHA(2), RHAT(10), RSL(NN), CHISQR 
B(M), KEY(M), KEY1(NN), KEY2(NN), XIJHAT(N,M) 

WEIB(M), BETNEW(N), XBENEW(N), XNEW(N), XIJNEW(N,M) 
LAMNEW(N), LAMBIG, RATNEW(N), LAMUNR, RSLNEW(NN) 
RHTNEW(10), CLR(10), CLRNEW(10) 

BHTBAR, BNRBAR, BHTMSE, BNRMSE 


S(N), DF, ICOUNT, SAMPLECN), CASE 


ISEED / 1736 / 
ALPHA / .2, .05 / 
RSYS / .90, .975 / 


BEM / jegri. &.106.408. 1.5, / 

BEM oe Wie Te Ope Steele oy Th ee SO oe Os 6d 
fey lent) 

| OO RI ORE) ORS) 

Dent Dis G DEG Ge: Oe Se OC mOr 2.8.9. 4. 
Dy 

x / 2,3,7,8,10/ 

2 S.7 8 0. 37.10, 1 em. LO 1n8y/ 


SAMPLE 
SAMPLE 
SAMPLE 


/ 97 ae. 6 Gy 
VS0R 21 A0n 567 
/ 30,63,75,98,26 / 


SAMPLE 
SAMPLE 
SAMPLE 


PO. 5 in Ord 0G oreo, 5 oe io Or 

3.7.10 1520S 67 Sa On On 10 47 lS ln e/ 
/15,40,35,17,26,67,50,65,80, 32,95, 100,15,45, 30/ 
SAMPLE / 15 * 15 / 


= 1 


Polo I = lt, N 


LE a SAMPLE( I) Pie TS) RoUBM( LT) = 1.40 
TEV CRSANE EE Gia he Je hoUb Cy «= "1-50 
CONTINUE 


pore Jon tea 


/* 
DO 


FINDING THE LAMBDA AND LAMBDA STAR FROM THE GIVEN DATA 
90 I=1,N 


* 
ve 


ris 


ve / 


ZU 


40 


60 


LAMBDA(I) = (( - ALOG(RSYS(JJ)) / N) ** (1. /BETA(I)))/ X(I) 
LAMSTR(I) = LAMBDA(I) ** BETA(I) 
sa = 0.84 * SAMPLE(I) - 0.5 

CONTINUE 


/* SUBROUTINE ZTABLE COMPUTES THE RIGHT PERCENT POINT ZALPHA */ 
/* FROM RIGHT CUMULATIVE PROBABILITY ALPHA o/ 


DO 40 I=1, 2 
CALL ZTABLECALPHAC I), ZALPHACI)) 
CONTINUE 


/* COMPUTE THE DEGREE OF FREEDOM FOR CHI-SQUARE ¥ / 


DF = 0 
pO 60 I=1,N 
DF = DF + SAMPLE(1) 
CONTINUE 
DF = 2 * DF 


/** COMPUTE CHI-SQUARE VALUE FOR GIVEN ALPHA AND DF ve / 
DO 200 i= ae 


IF ( DE . EGae] THEN 
P = ALPHACLL) / 2 
CALL ZTABLE(P, ZALPHACLL) ) 
CHISQR = ZALPHACLL) *** 2 


ELSE Bate) DF = i072 ees. 
CHISQR = -2 * ALOGCALPHACLL)) 


ELSE IF ( DF .GE. 3 ) THEN 
CHISQR = DF * (1 - 2./(9 * DF) + 
ZALPHA(LL) ** SQRT( 2./(9 * DF))) ** 3 
ENDIF 


LAMMAX 


0 
LAMBIG 0 


DO 180 ITER = 1, NN 
DO 140 I=1, N 


/* GENERATE EXPONENTIAL RANDOM NUMBER WITH MU=1 */ 
CALL LEXPNCISEED WEP Ooo Li ely ele ) 


/* CONVERT EXPONENTIAL RANDOM NUMBER WITH MU=1  */ 
/** TO MU=LAMBDA AND GET WEIBULL RANDOM NUMBERS */ 
DO 80 J = 1, SAMPLE(I) 


XX(I1,J) = EXPO(J) / LAMSTR(I) 
WO1.J) = XX(1.3) + Clearer) 
Y(I1,J) = ALOG(W(I,J)) 

WEIB(J) = W(I,J) 


80 


100 


130 


140 


BC J) = Y(I,J) 
CONTINUE 


/* SUBROUTINE SHSORT WILL SORT ARRAY B IN 
/* ASCENDING ORDER BY SHELL SORT ALGORITHM 
CALESHSORICER SeKEY, SAMELECID® 

/* FINDING BETA HAT (BETHAT) BY THE FORMULA 


SUM1 
SUM2 


HW ol 


0 
0 


DO O0e) = esc) 
SUM1 = SUM1 + B(J) 
CONTINUE 


DO 120 J = S(I)+1, SAMPLE(I) 
SUM2 = SUM2 + B(J) 
CONTINUE 


DECAL Ie =—eCOSMPLECI ys hSUBNC I) 7 
COC Ts SUMNZ 7 eC SAMPLEGI )S—es Gh j= sUNl) 


TEMP = BETA(I) 

/** SUBROUTINE NEWTON WILL COMPUTE THE ESTIMATE 
/* OF BETA (BETNEW) BY NEWTON-RAPHSON METHOD 
CALL NEWTON(WEIB, SAMPLE(I), TEMP, BETNEW(1)) 


XTOBET(I) = XCI) ** BETHAT(TI) 
XBENEW(I) = XC1) ** BETNEW(T) 


XHAT(1) 
XNEW( 1) 


oil 
Lae i SB) 


DO 130 J = 1, SAMPLE(I) 
XIJHAT(1,J) = WC1,J) ** BETHAT(I) 
XHAT( 1) XHAT(1) + XIJHAT(I,J) 


PA JINEWC 1 ae Wl SER TNEWCI) 


XNEW( I) XNEW(1) + XIJNEW(I,J) 
CONTINUE 
LAMHAT(I) = SAMPLE(I) / XHAT(I) 


LAMNEWC I) = SAMPLECI) / XNEWCI) 


IF (LAMHAT(I) .GT. LAMMAX ) THEN 
LAMMAX = LAMHAT(1) 
ENDIE 


IF (LAMNEW(I) .GT. LAMBIG ) THEN 
LAMBIG = LAMNEW(I) 
ENDIF 


CONTINUE 


= 0 


oy) 


7% / 
PAY 
ate if 
#% 


SUM3N 
SUM4 


0 
0 
SUM4N 0 


DO 160 I=1, N 
RATIO(I) AMHAT(I) / LAMMAX 
RATNEW(I) = LAMNEW(I) / LAMBIG 
SUM3 = SUM3 + RATIO(I) * XHAT(I) 


Ske 


SUM4 SUNG = RATIO LT) = TORE rar) 
SUM3N = SUM3N + RATNEW(I) * XNEW(1) 
SUM4N = SUM4N + RATNEW(I) * XBENEW(1) 
160 CONTINUE 
LAMMU CHISGR Sy (825 sul) 


RSLC ITER) = EXP( - LAMMU * SUM4 ) 


LAMUNR CHISOR>/ G22. Sule) 
RSLNEWC ITER) = EXP( - LAMUNR * SUM4N ) 


180 CONTINUE 


CALL SHSORT( RSL, KEY1, NN) 
CALL SHSORTCRSLNEW, KEY2, NN) 


KK = (1. - ALPHACEL) ie NN 


“KHAT IS THE KELTABIDITY OF THE SERIES SYSTEM COMPUT 

* BY THE FORMULA 

* RHTNEW IS THE RELIABILITY OF THE SERIES SYSTEM COMPUTED 
/* BY THE NEWTON - RAPHSON METHOD 


SS SS 


RHAT( CASE) 
RHTNEW( CASE) 


RSL( KK) 
RSLNEW( KK) 


/* SUBROUTINE -PINDI WILL SIND Tibet RODS CONT IDENGE Siler 
CALL FINDSCRSE NN, Ksvo Jape.) 
CALL FINDJCRSLNEW, NN, RSYS(JJ), JZ) 


CLR( CASE) 
CLRNEW( CASE) 


Ji / FLOAT(NN) 
J2 / FLOATC(NN) 


CASE = CASE + 1 


200 CONTINUE 


220 CONTINUE 
DO 240 I =1, N 
WRITE(6,500) BETA(I), SAMPLE(1) 
240 CONTINUE 


WRITE(6,550) N, NN 
WRITE( 6,600) 


38 


« af. a 
oe 


Se as a 


af 
Ped 


ve / 
ray 


ey eae a es ep 


WRITE( 6,650) 
WRITE(6,670) 


CASE = 1 
WOU ZOO aan, 2 
DO 260 LL = 1, 2 


WRITE(6,700) N,RSYS( JJ) ,ALPHA( LL) , RHAT( CASE) , RHTNEW( CASE), 
0 CLRCCASE), CLRNEW( CASE) 
CASE = CASE + 1 


260 CONTINUE 
280 CONTINUE 


500 FORMAT(T10,' INITIAL BETA grec) Pat alae 
vs Ras. NUMBER OF SAMPLE SIZES 2%, 13) 
550 FORMAT(T10,' NUMBER OF COMPONENTS :* 9X13, 
vs /,T10, "NUMBER OF REPLICATIONS a ls) 
600 FORMAT('1' 5(/), 7%, lest’ RELIABILITY OF SERIES SYSTEM sve" ) 
650 FORMATC/// 735, Rl, CLI - USING THE FORMULAR' 
ie TES, R2, CL2 ; USING THE NEWTON RAHPSON METHOD ' ) 
670 FORMAT(///_ T6, NUMBER OF' ,T19, ‘RELIABILITY’ ,133, ALPHA',T46, R1' 
ve T56, R2! ,T66, ' TRUE ,T74, 'TRUE' / 
vs T6, COMPONENTS ' HO Ol Svenel< Loom Gab i . Ty Gey Gaines: 
vs /,T5,76C' 
700 FORMAT(/,T8,15, 122, P5005 Pos,h oe, Ios ooo eo 


ee ore 


afealealontsa'satasr’s wlon'’s vlaata ~ wloaleatoatsaleatsniontsctantsatsats slesevlesleslesovedesesleslesesleste sie wloalsalanisatsatsatasat wl salealtsatonatsn!, alsatsaltontsntonlanisats alain a tn an a ate ea a 
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SUBROUTINE ZTABLE( ALPHA, ZALPHA ) 


<< SUBROUTINE ZTABL2 COMPUTES RIGHT PERCENT POINT ZALPHA FROM >> 
<< RIGHT CUMULATIVE PROBABILITY ALPHA >> 


REAL ALPHA, ZALPHA 


Wee CALPHA .GT. 0.0) .OR. CALPHA .LT. 1)) THEN 


W = - ALOG(4 * ALPHA * ( 1 - ALPHA )) 
ZALPHA = Ss aanae * (2.06118 - ( 5.72622 / (W + 11.6406)))) 


DE MCAEPHA .GT.. 0.5) THEN 


ZALPHA = - ZALPHA 
ENOL 


39 


ENDIF 


RETURN 
END 


Sededesesesccevesevesevedescdcolcskbvestcvevescseok teokacsesesdcslesdcsleclcdkvaesecdlesdcdesledesedkeskvskcaesleakleslevedesdvococvooracdedede dese deskesevesc se teste 


hy 


*, 


* 


Pr 


sesleseslesieslesisslevesecscolesisvesicdescceclesicseskatslestdeslcceslcotevestecesleceseskcdeakcedcscdesicsestsccvcslecleslcvlescckscccccsleseskestceskoledesssleskstcokeclese 
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SUBROUTINE SHSORT(A,KEY,N) 
DIMENSION A(N),KEY(N) 
M1=1 

M1=M1*2 

IF (M1 .LE. N) GO TO 6 
M1=M1/2-1 

MM=MAX0(M1/2,1) 

Gonoeed 

MM=MM /2 

IF (MM .LE. 0) GO TO 100 
K=N-MM 

DO 1 J=1,K 

II=J 

IM=II+MM 

IF (ACIM) .GE. ACI1)) sco 10m 
TEMP=A( 11) 

IT=KEY(II) 

A( I1)=AC IM) 
KEY( 11 )=KEY( IM) 

AC IM)=TEMP 

KEY( 1M)=IT 

T1=11-iM 

mel, (iy O) We) Ge 1 
CONTINUE 

GO TO 20 

RETURN 

END 


wean an celan's el snlsn’ wn se'salsninaionion's sieve eaten’ sales o's se aS, sesJe s‘es'esle wa ale utente aleulea'satan's ale els ute ale stevlenten! oe mle ala nts alee e's wisn sales nlentsalea'antanianien'cnwa'an's a), 
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SUBROUTINE NEWTON(TIME, M, BETA, BETNEW ) 


<< SUBROUTINE NEWTON WILL COMPUTE THE ESTIMATE OF BETA (BETHAT) >> 
<< BY NEWTON - RAPHSON METHOD Ba 


REAL TIMECM) , ENT IMB GROG) a ACOs 
I NIEGE Rees 


DATA TS / 6 / 


C = 0 
ITER = 


40 


R 
LNTS 


M 
ALOG(TS) 


DO. 20 T=" 1, R 
LNTIMEC I) = ALOGCTIME(I)) 
C = C + LNTIME(CT) 
20 CONTINUE 


C=C /R 


BODO 60 J= 1, 3 
SUM = 0 


DO 40 K=1, R 


IF ( J .EQ. 1 ) THEN 
SUM = SUM + TIME(K) ** BETA 
ELSE 
SUM = SUM + ((TIME(K) ** BETA)*(LNTIME(K) 2% (J-1))) 
ENDIF 
40 CONTINUE 


ACJ) = SUM 
60 CONTINUE 
foe FUNCTION PPRIME IS THE DERIVATIVES OF FUNCTION F ay 
OUOT = A(2) / AC 1) 
Pein = AC3) / AC 1) = CQUOT**2) 4+ ((1/BETA)**2) 
eeeoUOT = (1 / BETA) = C 


/* BETA IS UPDATED EACH TIME AND CHECK IF IT CONVERGES ve / 
BETA = BETA - F / FPRIME 


Pbk = ITER +. 1 
Pe CBETA .~GT.- 25.) GOTO 100 


fee ABS(F) .GI. 0.0001 ) GOTO 30 


ALPHA 
BETNEW 


(AC1) / R) ** (1/BETA) 
BETA 


RETURN 


100 WRITE(6,*) ‘DID NOT CONVERGE' 
WRITE(6,*) ‘TRY AGAIN WITH BETTER ESTIMATE OF BETA' 


RETURN 
END 


Neveu ew ON eu ON ON OU ANON ON OS GC OLEN GN ON EN ED EN AN ANON RODIN ERED RECON GN EN FTES FOES ON FO ENTE ICES ENV TN INTRON ED ON OD EVOL ON ON ONTO ELEN TE ON OO ON TL ON URES 
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SUBROUTINE FINDJ(CA, NN, R, J) 


<< SUBROUTINE FINDJ FINDS THE INDEX OF ARRAY A WHICH THE 
<< VALUE OF I] Is CLOSEST 10 5: 


REAL ACNN) , Ry Vanur 
INTEGER J 


VALUE = ABS(A(NN) - R) 
DOGO. = NNei eel, et 


IF (ABS(AC(I) - R) .LT. VALUE) THEN 
VALUE = ABS(ACI) - R) 


ELSE 
Nii ecm 
RETURN 
ENDIF 
CONTINUE 
RETURN 


END 


>> 
>> 
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APPENDIX D. FORTRAN CODE FOR INTERVAL ESTIMATION 


PROG 


PROCEDURE - NORMAL CASE 
RAM NORMAL 


devededevevedevededtededededcsedededededcdcacdescdedcdedededcvedcdedkvedcdesedvedcdevedededcdedetekdetedckkdedtedededededodeds vede 


ve 


THIS PROGRAM DETERMINE THE ACCURACY OF AN APPROXIMATE 
LOWER CONFIDENCE BOUND FOR PC X > Y ). 


PROGRAM NORMAL IS THE CASE WHEN Y IS GIVEN A VALUE YO 
& TA IS9NORHALLY DISTRIBUTED WITH 
UNKNOWN PARAMETERS. 


ssisdedodededcdkededededetedcdedkedededededtedesededcdedededededetoeddesevoddevededkededovededcdedcdvdededcdevekdededededssodededs deve 


REAL 
REAL 
REAL 
REAL 
REAL 
REAL 
REAL 
REAL 


INTE 
DATA 
DATA 
DATA 
DATA 
GAoE 
/* if 
DO 5 


/* 
ied 


ae 
/* 
CA 


fa 
#e 
is 

ar 


YO, P, ZP, MUX, SUMK, SUMX2, XBAR 

TEMP1, TEMP2, TEMP3, SIGHAT 

R(3), SIGMAX(3), ALPHA(2), X(100), Y(100) 

RL1(54,1000), ARL1(1000), BRL1(54), KEY(1000), CIRL1(54) 
TIRL1(54), T1C11(54), T2RL1(54), T2C11(54) 

TRIgn@o ee UGE 12) 

ZHAT( 1000), ZKNIFE(54), ZVAR(54), ZBAR(54) 

MUJUNK( 54) 


Gi NUiMxCs 5  NUMYC@3)), Chol... OUNE 


Ry .95.99...995/.. SIGNAGE) 1-207 Uns 10525.75/ 
NUMY /10,25,75/ , ALPHA / .2, .1/ 

ISEED, JSEED / 4875,7981 / 

YO / 400 / 
= 0 


[SS teeth OR Romma) =. 95, RO2) = -99, RUS) = .995 
OSS 


FP ISo0UHE RIGTHCUPPER) CUMULATIVE PROBABILITY 
= 1 SeRChL) 


SUBROUTINE ZTABL2 COMPUTES RIGHT PERCENT POINT ZP 
FROM RIGHT CUMULATIVE PROBABILITY P 


ibEeaUABL2( P, 2P ) 


Delo JHE INDEX TEOR SiGe Or x. 
SIGMAX( 1) = .5, SIGMAX(2)9=— 1, SIGHAKC3) = 3 


43 


ve / 


% / 


% / 
ve / 


ah al, 
~~ ~, 


C2 C2 


© 


50 


DO 400 JJ = 1, 3 
MUX = ZP * SIGMAX(JJ) + YO 
/*= KK IS THE INDEX OF NUMBER OF Xx. 


f* NUMRCT) = 10, NUMKCZ) = 25, NUNC] ee 
DO 300 KK = 1, 3 


/* LL IS THE INDEX FOR ALPHA. 
/* ALPHAC1) = .2, ALPHA(2) =.1 
WOO IS 1, 2 


CALL ZTABL2(ALPHA(LL), ZALPHA) 


CASE = CASE + 1 
TEMP = (NUMXCKK) = 1.) / NUMXCKkp 
TEMPO = SQRT( TEMP) 


MUJUNK( CASE) = MUX 
/* REPLICATE 1000 TIMES FOR BAGH CASES 
DO 100 I = 1, 1000 


/* USE NORMAL RANDOM NUMBER GENERATOR TO GET 
/* NUMXCKK) NUMBER OF xX. 


CALL LNORM( ISEED,X,NUMX( KK) ,2,0) 

/* THIS PART IS TO GET SAMPLE MEAN(XBAR) AND 
/** SAMPLE VARIANCE( XVAR) 

DO 50 MM = 1, NUMX(KK) 


XC MM) SIGMAS(C II )) 7 ACh) aa. 
CONTINUE 


/* SUBROUTINE VAR WILL COMPUTE SAMPLE MEAN AND 
je SAMPLE VARIANCE 
CALL VAR(X, NUMXCKK), XBAR, XVAR) 


/* NOW WE COMPUTE THE LOWER CONFIDENCE BOUND 


TEMP1 = (XBAR - YO) / SQRT(XVAR * TEMPO) 
TEMP2 = 1. / NUMX(KK) + 

((XBAR - YO) ** 2) / (2*(NUMX(KK)+1)* XVAR) 
TEMP1 = TEMP1 - ZALPHA * SQRT(TEMP2) * TEMP 


/* SUBROUTINE ZTABL1 COMPUTE RIGHT CUMULATIVE PROBA- 
/* BILITY FROM KiGHhI PERCE ie Oe 


CALL ZTABL1(TEMP1, ARL1(I)) 


sf Se 
™~, M 


7 / 


ARLIC1) = Tee ARLICT) 
CALL JKNIFECX, NUMXCKK), YO, ZHAT(I)) 
100 CONTINUE 


CALL VAR(ZHAT, 1000, ZBAR(CASE), ZVAR(CCASE) ) 


SIGHAT = 1. /NUMX(KK) + ZBAR(CASE)%*2 / (2%(NUMX(KK)+1)) 
SIGHAT = SQRT(SIGHAT) * TEMP 
SIGHAT = ZBAR(CASE) / TEMPO - SIGHAT * ZALPHA 


CALL ZTABLICSIGHAT, ZKNIFE( CASE) ) 

ZKNIFEC(CASE) = 1 - ZKNIFECCASE) 

/* NON-IMSL LIBRARY 'SHSORT' WILL SORT ARL1 BY SHELL ooh 
/* SORT ALGORITHM. uy 
CALL SHSORTCARL1, KEY, 1000) 

POmioO tt = 1, L000 


RLICCASE,I) = ARLICI) 
150 CONTINUE 


/* SUBROUTINE FINDJ FINDS THE INDEX OF ARL20 WHICH THE */ 


poe enue OF Ir LS CLOSEST T0-K. / 
GALE SEENDICARLI; RCI) 2J3) 

/* DIVIDING THE INDEX BY 1000 WILL GIVE US TRUE we / 
/* CON ERIDENCE LEVEL. ve / 


CiRbLCCASE) = J / 1000: 


MM = 1000 * (1 - ALPHACLL)) 
BRGICGASE) = RLI( CASE, MM) 


CALL TNYVALCR(II), SIGMAX( JJ), NUMXCKK), 


* ALEWAGEI) ,.. ZALPHAe mua, YO, TREI, TCI1) 
TIRLICCASE) = TRL1(1) 
Tila CeOASE) ==) Teli 1) 
M2Ze bc Ash)- = TRLIC2) 
Peer TCCASE) = Telic2) 
200 CONTINUE 
300 CONTINUE 


400 CONTINUE 
5900 CONTINUE 
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WRITE(6 ,600) 


WRilliCo,es 0) 
I= 1 
LINE = 0 


DONS > 02h 
DO° 350-410 2= se 
DOS SOK at, 3 


PE Ci NES oh DEN 


WRITE( 6, 800) 
LINE = 0 
ENDIF 


WRITE(6,710) R(I1), SIGMAX(JJ), NUMX(KK), 
vk BRLI(1), BREICI+1), CIRINC1), CIRmmCItt). 
vs ZKNIFE(1), ZKNIFE(I+1), 
MUJUNK(I), MUJUNK(I+1) 
eS ee 
LINE = LINE + 1 
550 CONTINUE 


WRITE( 6,610) 
WRITE(6,650) 
T=1 
LINE = 0 
DO 56011. — mies 
DO 560 JJ = 1, 3 
DO SoOuKk = aes 


TP GOUNE Gres. een 
WRITE( 6,800) 
LINE = 0 

ENDIF 


WRITE(6,700) R(II), SIGMAX(JJ), NUMX(KK), 
2 TIRLI(1I), TiRbI¢@+1), TICIiG@).. Tien ic ia 
IT=I1+2 
LINE = LINE + 1 
560 CONTINUE 


WRITE( 6,620) 
WRITE(6,650) 
I=1 
LINE = 0 
DO 570 11) =e 
DO) 5700 Jae 
DO S709 = Ios 
IF (LINE .GE. 3) THEN 
WRITE(6, 800) 
LINE = 0 
ENDIF 
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(Se al eG | a 


WRITE( 6,700) (RG jes OMAN Open Chk ) 
cl TRG. T2R pice eee Ol LCL). PZELIC I rl) 
tea 2 
LINE = LINE + 1 
p77 OS CONTINUE 


600 FORMAT('1',5(/),7X, '*** Y IS GIVEN A VALUE OF YO *#%r') 


og) FORMAT( i ,5(/),7X, ‘deere =Y IS GIVEN A VALUE OF YO **e*' , 
[/ ,20X,' -- X IS TRUNCATED NORMAL WITH PROBABILITY 5a 5 

620 FORMAT’ 1 »5¢/), 7X, Nieves = Y TS GIVEN A VALUE OF YO 9 *tev! 

//,20X, X Is TRUNCATED NORMAL WITH PROBABILITY .90 --') 
630 FORMAT(3(/),5x, CASE ',6X,' RL1(1000*(1-ALPHA) ) i 

* 4X, 'TRUE CONFIDENCE LEVEL’ ,5X,'JACKNIFE' , 

SOR SX N',3X,' ALPHA = .2' ,12X, a ke 2 

* Wee at Oko = ONG 5 = ak, BC So ox rie ae vO 

650 FORNAT(3(/) Ep.g CASE ,6X,° RL1(1000*( 1-ALPHA) ) me 


4X, TRUE CONFIDENCE LEVEL ik a, *R SX N’ 3K, 


"ALPHA = .2' Wey LU Sar CeCe ex 
PG ye ake 1C 7) 


700 FORMAT(3X,F4. 3,2X,F3. 1,2X,12,4(7X,F7.4)) 
710 FORMAT(3X,F4.3,1X,F4.1,2X,12,4(7X,F7.4),1X,2(2X,F6. 4) ,2(2X,F9.5)) 
800 FORMAT(/) 


RETURN 
END 


we elon on sole nt selentanlsn'we'eusalanlsalsdewselsuiseswseseisetsstsulententsuiscisoleriscisnltentsaeatsedevsatedodsetocisneuteaieesalee'suenu'esnientac’sutentantan'satnan'en sutan'scton ante 
VOTE GES EOS GE GN GN ON aN oN Oe AN OR Ae Ce EE RES CO ERODED TESS EV ODEN ED ED FS EDEN ED EDEN EN FR EN OR ED FN EN GS ENED ENTS EDO EN CL EN ES ED ENED ENED ON ODEN OS 


. ale sie eos we we we ae ate aloals Se ale ada we ale ale ule ale ale aleu's we wa ate we we ale ale a's we ale we ele aloule ale afe ae ale ale we ale ule afe ae ae ae we ale ale we ate ala ale AU eta we oe ate ae 
Cu GS EV EDV ED ED AY CX FH EH FLFR ER FR FVD ED EN GH ER GV GL GV ERM EX EX FL TDL GX FS £9 FR FX EY EX AR GH GV ERGY GX EX FX EX FR EX EYL FX FH FH FH EX ER Gd FX EH CRO FR SY FH 


SUBROUTINE ZTABLI(ZALPHA, ALPHA) 


<< SUBROUTINE ZTABL1 COMPUTE RIGHT CUMULATIVE PROBABILITY FROM >> 
Seek IGhl PERCENT POINT. >> 


REAL PI, ALPHA, ZALPHA 


PARAMETER {( PI = 3.141592 ) 


ie ZAGPHA .EQ. 0.0 ) THEN 


ALPHA = 0.5 

ELSE 
ALPHA = EXP( -2 * ZALPHA ** 2 / PI ) 
ALPHA = SQRT(1 - ALPHA * (1 + 2*(PI-3)*ZALPHA**4 /(3*PI%*2))) 
edee="0.50-°( f | ALPHA ) 


eC 4AceHA . bf — 0.0 ) THEN 
ALPHA = 1 - ALPHA 
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BND LE 
END 


RETURN 
END 


tof Tantlactaa! octactsciact a selanlactantactae! 1 Weta! Ae $3 cto ces a 
Sevedevetevevesevevevoevesevetesevevesedededetedesvesedededsededsededesedetetesedededekdedesedededesededsedesevedesedededededesedesededededesesesede 
ate Sa 
ra ae 


UFC EL EN EU GU EU EL FR EULER ER FL FY od FB ee GB GP Fe CY FB FU FT. Ed HR EH FN 7 vs ev ev. ar ee wes x sR ed ev FR wee er varks PR ER EL EDT EL CD EN 4 varxy Pieris edu ds rAerty ey ee rieris eves 


SUBROUTINE ZTABL2( ALPHA, ZALPHA) 


<< SUBROUTINE ZIABL2 COMPUTES RIGHT PERCENT POINT ZALPHA FROM ee. 
<< RIGHT CUMULATIVE PROBABILITY ALPHA >> 


Sp Sp a ae xe 


REAL ALPHA, ZALPHA 
IF ((ALPHA .GT. 0.0) .OR. (ALPHA .LT. 1)) THEN 


W = - ALOG(4 * ALPHA * ( 1 - ALPHA )) 
ZALPHA = SQRT(W * (2.06118 - ( 5.72622 / (W + 11.6406)))) 


IF CALPHAS. GT, US ieethEN 
ZALPHA = - ZALPHA 
ENDIF 
EINES 


RETURN 
END 


SRN GN EU ER EL EOC ER IN GR FUGUE TROD EN GO GN ED EC ES RELIED EVES ECAR ECON EN TREC PREV TREN EL ERT EN GN ODER ES FLEE FL INGTON OLIN FUER FLINGER TR ELEVEN CO EL IN ED ES ON GN OD 
ale 
At 


SU GR EV EU ER FEUER ED GEL FR TE FR GE FR FU FU ER EVID FUL GB EL FRR EN FR FR FR FR GU FU INU GU ED EYL IH FR FU FR GL GL FU FR FUER EB EY FN GY ER FRR IH FR EL FER ER ITED ER FR ER IDV GVTEV EV ER GR ER Ed FN 


SUBROUTINE TTABLE( ALPHA, NU, TALPHA) 


<< SUBROUTINE TTABLE COMPUTES RIGHT CUMULATIVE PROBABILITY >> 
<< TALPHA FROM T-DISTRIBUTION BY GIVEN ALPHA AND NU >> 


CS a A eG 


REAL ALPHA, ZALPHA, TALPHA, Al, A2, A3 
INTEGER NU 


IF ((ALPHA .GT. 0.0) .OR. (ALPHA .LT. 1)) THEN 
CALL ZTABL2(ALPHA, ZALPHA) 
CANTER EN ee 


(S*ZALPHAT*4 + 16*ZALPHA**2 + 3) / 96 
(3°*ZALPHA***6 + 19*ZALPHA*4 + 17*ZALPHA**2 - 15) / 384 


a 
No 
Hon TT 


TALPHA = ZALPHA * (1 + A1/NU + A2/NU**2 + A3/NU**3) 
IF € (APHAUe Gls) 005 9) HEN 


TALPHA = - TALPHA 
ENDIF 
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ENDIF 


RETURN 
END 


Secvededvdevevevesedevededesesedededededcdvakvvdededededededcdcdedesesededededesevedededcdkvveskcsededetedcdededcdvdededovededeskdedtdevededede sede 
ve se 


Seavseslesevededestevedevedeslcdestevesteseveakvovlevesleslevlcscslevleseslestesdleskeseslcaksicseulcckakvesieveseakseslesesesledcckakesesevestdedlevestevlesleste clevieveves'c 


SUBROUTINE SHSORT(A,KEY,N) 
DIMENSION ACN) ,KEY(N) 
M1=1 
6 M1=M1*2 
IF (M1 .LE. N) GO TO 6 
M1=M1/2-1 
MM=MAXO(M1/2,1) 
eer 16 21 
20 MM=MM/2 
IF (MM .LE. 0) GO TO 100 
21 K=N-MM 
peeO 1 J=1,K 
ll=J 
11 IM=II+MM 
WoMCoCiy) GE. ACII)) GO TO 1 
TEMP=A(I1) 
FT=KEY(11) 
AC II)=A( IM) 
KEY( 11 )=KEY( IM) 
AC IM)=TEMP 
KEY( IM)=IT 
ITI=II-MM 
fe tt. cT. 0) GO TO 11 
1 CONTINUE 
ee.10 20 
100 RETURN 
END 


wsalsatoats a'snlen'on'e wlousalealtsatasulsatsctsalonten'esutna'sasatse'sulselsuisaiselseseioatsclisutselsuts ole ulse'su'sviselseisatsciselssivu'selsntselssisa'setselsalonisealsatscianlsulents als ulsn'selsates afoats wloute 
FR ER FR ERED ER FR ER ER IR EL FR EV ED FL ELEN FR ERED ER ER FR FL ER ER GR EL ER FH FR ER IR ES FR IL ERIS FH EL TL EH ED IN FD FD Aerny rer hs FV EV EV EL ER EDL FRB ER IL ER EVR EN FBR ER ER TH Hdd sed “% 


Saeiselsnentsnteatse'entseisalsetseoviantsaeelscisvian'snievissians a'selawlanlsatevisulouts wlaaltan's afealsatea’setsuisalsnten'salselsn'sulsulsalsn'setsa’se'se'sr sulavisalsalsn'seise'se'selse'se'selsn’onlon'ou'es 
NR eR ese er. ed. ev. er fu ert @ VER ER EDL ER ER ER EVER EV EV ER ER EN ER EL ER EN FR ER PR ES FS FRED GER odo GR ER HB E> FR ER ER EN Fd FD rierks REL EL EN EN FN CD ER ER EV EX EX EX FX EX ER IN CRED IN FR EN ED 


SUBROUTINE VAR(Z, NUMZ, ZBAR, ZVAR) 


<< SUBROUTINE VAR WILL COMPUTE THE MEAN AND THE VARIANCE >> 
a= OF ARRAY Z >> 


Sp Te Gp i Gp | 


REAL Z(NUMZ), SUMZ, SUMZ2, ZBAR, ZVAR 


SUMZ 
SUMZ2 


0 
0 


DO sO Our: 


1, NUMZ 


SUMZ 


SUZ eee ol) 
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CoC) C2 C2 


artealealeaanalselsnisatsetaasaisatsntselealeuise'saiasntsalsntsats lone ata ates o's ates ealasale atentoelsalselsaisaisalselaeisalselaniaets ulsalsa's afante alsaleatsatsalsn 
@e eR ee ey es 4d #d RV EV AR 4H FR FY. FY FRX EX GX #Y #9 FR Fd AN AN 42% FH FF riety OV EN EV EN ER AY FV FY ES GX GX 4B ER GOAN REV EN ER ER GV TN EX HY EN 


C2.C2 C2 <9 


SUMAZe= SUNZZ eZ 1% Z 
100 CONTINUE 


ZBAR 
ZVAR 


SUMZ / NUMZ 


RETURN 
END 


feutsa!, 


SUBROUTINE FINDJCA, R, J) 


ABS( SUMZ2 / NUMZ - ZBAR ** 2 ) 


<< SUBROUTINE FINDJ FINDS THE INDEX OF ARRAY A WHICH THE 


<< VALUE OF I7 15 CLOSES? 103) 


REAL A( 1000), R, VALUE 
INTEGER J 


VALUE = ABS(A( 1000) - R) 
DO 00s Seo lie 


TF CABS( AGI) ~- 3) = Li. VALUE = THEN 
VALUE = ABS(ACI) - R) 
ELSE 
Jean el 
RETURN 
ENDIF 


100 CONTINUE 


RET” RN 
END 


SUBROUTINE TNYVAL(R,SIGMAX , NUMX , ALPHA, ZALPHA ,MUX, YO, 


os AERC. Tere 1 ) 


<< SUBROUTINE TNYVAL IS THE CASE WHEN X IS TRUNCATED 


<< AND Y IS GIVEN A VALUE 


REAL R, SIGMAX, ALPHA, TALPHA, MUX, P(2), ZP(2), A(2) 


KEAL, AGICG) a, EP ee eit 


Pars) 


alas es 


ate ale ate ale aloala ole aloalante 
RP ev ery BEX EX FRY FO Fd 4% 


NORMAL 


REAL TARLIC1@G0), KEYC10G0) > TRLIC 1000) TBRETGZ Ree TCiRT Gs 


INTEGER NUMX, CASE, COUNT 


La a) 
>) 


ale le 
¢ 


estesesteslesiesievevevevevesesedevededesodestededededtedededededesesedededesedesedtededsedededededededededesededededolesdededkdedetedtedesesvedtedededsede te 


Je 


GRR EN CL ON CV ON EN COED EV TED FR FLERE SO BL EN ERED EC ERED EN ERED EV EDEN ED ER ES ENGR ED ELEVEN EN ER ELEN ER ENERGON EER ER EN OC OVEN ER EN ENON EN EN SO ER EN EN ER OREN GR ED 


>> 
>> 


alee wees als alan salen ante alsalontonts alaulaalaalaatsatselentsats uleusale alealsalsaleats Jeeves ale ale a alan'’s alsalon's vesles'e wlanleclsaalsciselseisctee'sn’s elses elses ale els u'se'selsaleuis ulveleeeatsn's 
ri eriy FR AD RR EL 49 YR OR FR #e OR AN ES 2H FR ES FR re Ad 4s 4e 8 4 riety Cv. egu deve > @% #0 €B #B 4% #8 #2 Fd ve #8 ah #0 FO FR FY FD FR FL £2 FR #9 FX FH FR 2H 49 FR FR FR FX FX FR riveris @#y #% 


ales 
@eN 


ate 
¢e 


>> 


100 


200 


DATA P / .95, .90 / 


WetAw tobe le JoBED / 4875 ..7961 / 


DO 400 I =1, 2 


Poh jo— =. PCT) 
COUNT = 0 


CALL, ZTABE2ZCP(1).,. Z2PCE)) 


AC1) = "MUX  ZP( 1) * SIGMAX 


DO 200 J = 1, 1000 
CALL LNORM(ISEED, X1, 


PeEMP SS .MUX + SIGMAX * 


Lee 2S Os) 


X1 


IF (XTEMP .LE. ACI)) THEN 


COUNT = COUNT + 1 
X( COUNT) = XTEMP 
ELSE 
GOTO 100 
ENDIF 


PESCCOUNT 2L1. NUMX) GOTO siG6 


CALL VAR(X, NUMX, XBAR, XVAR) 


/* NOW WE COMPUTE THE 


CEE 
TENP2 = 1.7 NUMX > 


(CABAR = YO) +" 


LOWER CONFIDENCE BOUND 


(XBAR - YO) / SQRT(XVAR*(NUMX-1) / NUMX) 


2) 7) (2 5CNUMNAF1) = KVAR) 


TEMP i= PEP L =—ZAGPnAY “SORT TEMP 2 ) 


/* SUBROUTINE ZTABL1 COMPUTE RIGHT CUMULATIVE PROBA- 


J BILITY PRC RIGHT, PERCENT £OInNt. 


CALL ZTABLI(TEMP1, TARL1I(J)) 


Tania I j= a. = PARED 


CONTINUE 


/* NON-IMSL LIBRARY 'SHSORT' WILL SORT TARL1 BY SHELL 


/* SORT ALGORITHM. 


CALL SHSORT(TARL1, KEY, 1000) 


DO 300 J = 1, 1000 
TRL1(J) = TARL1(J) 


a 


= 


300 CONTINUE 
G /** SUBROUTINE FINDJ FINDS THE INDEX OF TARL1 WHICH THE 
C /* VALUE OF IT IS CLOSEST TO R. 
CALL FINDJ(TARL1, R, J) 
C /** DIVIDING THE INDEX BY 1000 WILL GIVE US TRUE 
C /** CONFIDENCE LEVEL. 
TCIRI(I) = J / 1000. 


MM = 1000 * (1 - ALPHA) 
TBRLIC1I) = TRL1CMM) 


400 CONTINUE 


RETURN 
END 


a 
Ts SS 


COVES EN GV EV EN CL ed AN eRe SN ON OV EN FRED ER ER EN GN CR Ee EDEN GLEN EO EN EO EN EDEN EV EVER ERED LEC EDEN ENON ELEN EL EV CO ELEN ELEN EDEN ES EN EN EN ES ELEN EN EN EVEN EL ES EVEN GS 


ver 


alentsefsaioe'saton'snlon'sulenlon'sn'sn'en'eulsulo se ale alantantanlasarlevedenlesliantenen'sclantsalentsnieviocoleuaslentonlooteniontscisntantantententenionselan's alsa'scaulodau'sate 
x. aerne #yx ey ¢% Berner ere ers #9 er eh #8 4 OR EV EV EVER EL EV ER FEVW ER CB FR EL FUL ER FR FR ER ER ER EYL EL FR FL EB ER EU’ ER 4% FH FR ER EV EX EY PL ER ER ER AX GX FN Cd 


SUBROUTINE JKNIFECK, NiO eZine) 


<< THIS ROUTINE IS FOR THE ‘JACKNIFE' METHOD. 
<< DUMMY PARAMETER IS ZHAT. 


2 Il Se Bal aT 


elouteulousutonten'sn'en'sn's 
RGR EV EX REN CL ER EV Gd 


>> 
>> 


REAL X(N) ,AOMIT( 100) , XOMBAR( 100) , XOMVAR( 100) , ZOMHAT( 100) , ZHAT 


ENG EGER bie iN 


pO 200 I = 1, N 
DO 100 J=1, M 


IF (J .GE. I) THEN 
XOMIT(J) = X(J+1) 
ELSE 
OOO I) = IH) 
ENDIF 
100 CONTINUE 
CALL VAR(XOMIT, M, XOMBAR(I), XOMVAR(I)) 


ZOMHAT(I) = ( XOMBARCI) - YO ) / SQRTCXOMVARCI)) 
SUM = SUM +ZOMHAT(T) 


200 CONTINUE 


ZHAT = SUM / N 


A” 
to 


RETURN 
END 
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